function [fars, frrs] = slhistroc(ha, hr, op)
%SLHISTROC Computes the ROC curve from value histogram
%
% [ Syntax ]
%   - [fars, frrs] = slhistroc(ha, hr, op)
%
% [ Arguments ]
%   - ha:           the histogram of the values that should be accepted
%   - hr:           the histogram of the values that should be rejected
%   - op:           the option of for thresholding
%                   \{:
%                      - high:  accepts when v >= threshold
%                      - low:   accepts when v < threshold
%                   \:}
%
%   - fars:         the false accept rates at the histogram bin edges
%   - frrs:         the false reject rates at the histogram bin edges
%
% [ Description ]
%   - [fars, frrs] = slhistroc(ha, hr, op) computes the ROC and
%     returns the false accept rates and false reject rates corresponding
%     to the thresholds set at the edges of histogram bins.
%
%     In the input, ha, hr are two histograms computed from the the values
%     that should be accepted and that should be rejected respectively.
%     Suppose va, and vr are lists of these values, then
%       $ ha = histc(va, edges);  hr = histc(vr, edges); $ 
%
%     In the output, the sizes of fars and frrs equal the sizes of ha and
%     hr.
%
% [ History ]
%   - Created by Dahua Lin, on Aug 8th, 2006
%   - Modified by Dahua Lin, on Jul 19, 2007
%       - make the semantic conform to histc
%       - refactor the program structure
%

%% parse and verify arguments

error(nargchk(3, 3, nargin));

assert(isnumeric(ha) && isvector(ha), 'sltoolbox:slhistroc:invalidarg', ...
    'ha shoule be a numeric vector.');
assert(isnumeric(hr) && isvector(hr), 'sltoolbox:slhistroc:invalidarg', ...
    'hr should be a numeric vector.');
assert(ischar(op) && (strcmpi(op, 'high') || strcmpi(op, 'low')), ...
    'op should be either ''high'' or ''low''.');
accept_high = strcmpi(op, 'high');

assert(isequal(size(ha), size(hr)), 'sltoolbox:slhistroc:invalidarg', ...
    'ha, hr, and edges should have the same size.');

if size(ha, 2) > 1
    original_row = true;
else
    original_row = false;
end


%% Compute ROC

if accept_high   
    frs = count_low_end(ha, original_row);
    fas = count_high_end(hr, original_row);
else
    frs = count_high_end(ha, original_row);
    fas = count_low_end(hr, original_row);
end

sa = sum(ha);
sr = sum(hr);

if sa > 0
    frrs = frs / sa;
else
    frrs = zeros(size(frs));
end

if sr > 0
    fars = fas / sr;
else
    fars = zeros(size(fas));
end


function x = count_high_end(h, original_row)

cs = cumsum(h);

if original_row
    x = [cs(end), cs(end) - cs(1:end-1)];
else
    x = [cs(end); cs(end) - cs(1:end-1)];
end

function x = count_low_end(h, original_row)

cs = cumsum(h);

if original_row
    x = [0, cs(1:end-1)];
else
    x = [0; cs(1:end-1)];
end

